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Localization behavior of vibrational modes in granular packings 
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Abstract. - We study the localization of vibrational modes of frictionless granular media. We 
introduce a new method, motivated by earlier work on non-Hermitian quantum problems, which 
works well both in the localized regime where the localization length £ is much less than the linear 
size L and in the regime f > L when modes are extended throughout our finite system. Our 
very lowest frequency modes show "quasi-localized" resonances away from the jamming point; 
the spatial extent of these regions increases as the jamming point is approached, as expected 
theoretically. Throughout the remaining frequency range, our data show no signature of the 
nearness of the jamming point and collapse well when properly rescaled with the system size. 
Using Random Matrix Theory we derive the scaling relation £ ~ L d ^ 2 for the regime £ 3> L in d 
dimensions. 
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Introduction. — Over the past years, many questions 
concerning the behavior of disordered systems have been 
put in a new perspective by addressing them from the 
point of view of the more general jamming scenario [1]. 
Especially for granular systems it has turned out to be 
very fruitful to study the changes in the properties and 
the response of granular packings as one approaches the 
jamming point from the jammed side, where the packing 
gets close to an isostatic solid. An isostatic packing is 
indeed essentially a marginal solid which has just enough 
contacts to maintain a stable packing. From simple count- 
ing arguments, one finds that the average coordination 
number Z of an d-dimensional isostatic packing of fric- 
tionless spheres equals Z lso = 2d [2]. Upon approaching 
this marginal solid, many static and dynamic properties 
exhibit anomalous behavior, associated with the fact that 
the excess number of average bonds, AZ = Z — Z iso , goes 
to zero [3-6]. In fact, AZ itself scales anomalously, namely 
as the square root of the difference in density from the one 
at jamming [3]. Likewise, the ratio G/K of the shear mod- 
ulus G over the compression modulus K is found to scale 
as AZ, and the density of states of the vibrational modes 
becomes fiat at low frequencies above some crossover fre- 
quency lo* ~ AZ, due to the emergence of many low fre- 
quency modes. Much of this behavior was explained by 
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Wyart et al. [4-6] in terms of the existence of an impor- 
tant cross-over length scale I* ~ 1/ AZ, the length up to 
which the response is close to that of an isostatic pack- 
ing. This scale £* diverges as the jamming point is ap- 
proached, but is difficult to probe directly. Nevertheless, 
the length £* has recently been uncovered as the impor- 
tant cross-over length to continuum behavior in the static 
response [7,8]. Although most of these results pertain 
explicitly to packings of frictionless spheres, there are sev- 
eral indications [9,10] that many of these observations and 
ideas can be generalized to frictional packings. 

It has been noted in several studies that both the re- 
sponse to a local or global deformation [7, 11] and the 
behavior of the vibrational eigenmodes [4, 6] of a packing 
become much more disordered as one approaches the jam- 
ming point : as the snapshots of two vibrational modes in 
Fig. 1 illustrate, far above the jamming point the eigen- 
modes have a structure reminiscent of what one gets in a 
continuum theory of an elastic medium, but close to the 
jamming point one is immediately struck by the appear- 
ance of many disordered "swirls" . The arguments put for- 
ward by Wyart et al. [4-6] indicate that the excess low fre- 
quency modes cannot be localized on scales < I* since they 
are the vestiges of the global floppy modes that emerge at 
the isostatic point. Hence, if there are any low-frequency 
modes away from jamming and if indeed their localization 
length is > £*, we should see this as the jamming point 
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Fig. 1: Snapshots of two low- frequency eigenmodes in our pack- 
ings. The arrows indicate the direction and magnitude of the 
displacements of the individual particles, (a) At low pressure 
p = 10 -6 , close to the jamming point, the mode is very disor- 
dered, whereas at high pressure (b) p = 3 • 10~ 2 , the mode is 
more reminiscent of an elastic shear wave. Similar features are 
seen in the response to a local or global deformation [4,7,11]. 

is approached. The aim of this paper is to investigate 
whether this is indeed the case. 

Localization was discovered fifty years ago by Ander- 
son [12], who in his study of non- interacting electrons in 
a random potential found that disorder can induce elec- 
tron localization. Unlike the extended (delocalized) Bloch 
waves, in a localized state the weight of the electron wave 
function is concentrated near some point in space; the am- 
plitude falls off as e~ r ^ with distance r from the center. 
This defines the localization length £,(E) which depends on 
the electron energy E. The possibility that disorder can 
localize the eigenmodes of systems governed by wave equa- 
tions is quite general and extends to many systems, not 
only sound modes [13-15] but also gravity waves [15], light 
propagation [15] and diffusion on random lattices [14,15]. 
We will focus on the localization behavior of vibrational 
modes of 2d frictionless packings. In two dimensions there 
is no localization-delocalization transition: in the presence 
of disorder the states are generally localized in the ther- 
modynamic limit for cases like the one presented here [13]. 

The dynamic response of granular packings is affected 
by three types of disorder — bond disorder, mass disorder 
and topological packing disorder. Any of these is suffi- 
cient to cause localization, but in practice all three play a 
role for realistic models of granular packings: bond disor- 
der is present for all force laws except one-sided harmonic 
springs, polydispcrsc particles will have varying masses, 
and topological disorder is naturally present except for es- 
pecially prepared regular piles, like a regular stack of mar- 
bles. Of course, in computer models these effects can be 
separated easily; we will not attempt to disentangle these 
three contributions here, but do use this freedom later to 
our advantage in testing our scaling predictions. 

The crucial dilemma in extracting the localization 
length of the vibrational modes of granular packings is 
that the effective disorder is so weak that one needs 
prohibitively large systems to reach the true localization 
regime £ <C L for most modes. Here L is the linear system 
size. At the same time, existing methods which are based 
on spatial averages (like the direct expression based on the 




Fig. 2: Scatter plot of the angularly averaged £'s of all the 2000 
modes of our granular packing of 1000 particles as a function 
of the frequency uj at pressure p — 4 • 10~ 6 studied with the 
method explained in the text. Note the large scatter and the 
fact that the £ values are of order of the linear system size 
L — 45 or larger throughout most of the frequency range. 

second moment of the eigenmode or the (Inverse) Partic- 
ipation Ratio method [16]) do not give much insight into 
the structure of the modes when £ approaches the system 
size L, i.e., for modes which are extended throughout the 
finite system. As Fig. 2 illustrates, this is the relevant 
regime throughout most of the frequency range, as only 
the modes with the highest frequency u> are truly localized. 
The method we introduce in this paper, which is motivated 
by earlier work on non-Hermitian quantum problems [17], 
is based on studying the response to an asymmetric per- 
turbation. It not only gives the proper localization length 
£ of each localized mode, but at the same time assigns 
a well-defined and precise direction-dependent value £(</>) 
to each mode, that spans through our finite system - 
see Fig. 3. We stress that although wc will follow common 
practice in referring to £ as the localization length even for 
£ > L, one should keep in mind that many modes extend 
throughout our finite periodic system, as both Figs. 2 and 
3 illustrate. As we shall show, this method does allow us to 
study the scaling with system size and disorder, and opens 
up the possibility to bring Random Matrix Theory [16, 18] 
to bear on this class of problems. 

While all methods essentially yield the same localization 
length in the localization regime £ <C L, the extension of 
the concept of a localization length to the regime £ ~ L 
depends on the method used and it is not a priori clear 
what £ in this regime pertains to. For our method, one can 
however extract useful information about the large system 
limit from studying the regime £ <~ L. In conventional 
methods, one finds £ w L if the system size is too small. 
With our method wc find a disorder-dependence too which 
can be used to extract quantitative estimates of the intrin- 
sic localization length. As we will discuss in a forthcoming 
paper, this is simplest in one dimension where we predict 
and find a scaling £ ~ AL 1 ! 2 in the regime £ *~ L. Since 
we expect a crossover to the localization regime when the 
intrinsic localization length obeys £j„t ~ L, this to es- 
timate the infinite size localization length from the small 
system data simply as £i„t ~ A 2 . Preliminary analysis [19] 
indicates that this simple estimate works well in Id, but 
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Fig. 3: Polar plot of the localization length in one of our 
granular packings at p — 4 ■ 10~ 4 at low in (a), high in (b), and 
intermediate frequency in (c). The angular variation of £(0) is 
comparable to the angularly averaged value itself. 

we focus here on the behavior as function of the distance 
from the jamming point in 2d. 

Our main results can be summarized as follows (i) The 
average localization length of granular packings is 

largely independent of the pressure, and hence of the dis- 
tance from the jamming point, (ii) However, away from 
the jamming point there are a few "quasi-localized" low- 
frequency modes which disappear when approaching jam- 
ming. This behavior is qualitatively in accord with the- 
oretical expectations for the change in behavior near the 
jamming point, (in) In accord with what is expected on 
the basis of Random Matrix Theory (RMT) [16,18], modes 
with £ < L are effectively noninteracting and the distri- 
bution of their level spacing is Poissonian, while modes 
with £ > L show level repulsion: the level spacing follows 
the so-called Wigner surmise of RMT. (iv) In the regime 
£ > L, £(u)) scales as L d / 2 and is inversely proportional to 
the disorder strength, in d dimensions, (v) Due to level 
repulsion the distribution P(£) falls off for large £ as l/£ 3 . 

Method. — We use 2d packings of 1000 frictionless 
particles which are prepared using molecular dynamics 
simulations — see [7, 9, 10] for the description of our al- 
gorithm that gently prepares packing at a target pres- 
sure and other details. The particles interact with the 
3d Hertzian force law, /y ^ 5y , where 6%j is the overlap 
between particles i and j. The unit of length is the av- 
erage particle diameter. Unless noted otherwise we here 
present results for our most extensive studies with 20% 
polydispersity in the radii, but runs with different amount 
of polydispersity give similar results. The masses of 
the grains are taken proportional to Rf , corresponding to 
packing of spheres in 2d. The confining pressure, with 
which we tune the distance from the jamming point, is 
in the range p G (10~ 6 ,3 ■ 10~ 2 ) in the units of the 
Young modulus of the particles. We employ periodic 
boundary conditions in both directions. Our use of the 
3d Hertzian force law implies that the vibrational bonds 
hj = dftj/dSij ~ 8\( 2 ~ p 1 / 3 are disordered (they vary 
from bond to bond) and get weaker at smaller pressures. 
The natural frequency scale therefore goes down with pres- 
sure as p 1 / 6 . As in [9], when reporting our data we will 
therefore always rescale all frequencies lo with a factor 
p -1 / 6 , as to be able to compare data at different p. 

The vibrational modes and their density of states (DOS) 
are obtained in the standard way, by expanding the en- 
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Fig. 4: (a) Histogram of the ratio of squared amplitudes of 
the fourth (quadrupole) and second (dipole) harmonic at p = 
4-10 -4 . Most modes have predominantly dipole symmetry, (b) 
Average angular anisotropy A£ /£ as a function of frequency for 
various pressures. 

crgy about the equilibrium positions of the grains up to 
quadratic terms. Just as in solid state physics, the dynam- 
ical matrix, whose elements are the second derivatives of 
the energy with respect to the positions of the grains, de- 
termines the linear equations of motion of the vibrational 
modes. The dynamical matrix of a granular packing is a 
sparse symmetric matrix, because each particle only inter- 
acts with a few others. 

Our method to extract the localization length is mo- 
tivated by the work of Hatano and Nelson [17] on the 
derealization transition in non-Hermitian transfer ma- 
trix problems arising in the statistical mechanics of vortex 
lines in superconductors. Consider first the case of a one- 
dimensional chain of masses connected by springs with 
spring constants kij (j = i =fc 1) and periodic boundary 
conditions. We introduce an asymmetric bias term into 
the equations of motion so that the eigenvalue equation of 
a mode Uie~ tult becomes 

m l uj 2 u l = ^2 k v (e hi ' SlJ Uj - m) . (1) 
j=i±i 

Here the rest positions of the particles and Sij is 

a vector pointing from particle i to particle j. For h = 
this is simply the dynamical equation for vibrations. The 
trick now is that we can extract the localization length 
of each mode k by following whether or not its eigenvalue 
u)\ changes when we turn on h in small steps. Indeed, 
as long as h < l/£& the eigenvalue u\ will not change at 
all. To see this, note that in this case we can perform 
a "gauge transformation" to a field Ui = uie hxi which 
obeys the original equation with h = and which falls 
off exponentially on both sides so that, in a large enough 
system, it obeys the periodic boundary conditions. This 
implies that for h < l/£fc, the eigenvalue uo 2 , does not 
change. However, once h > ^ the function u obtained 
with this transformation does not fall off exponentially to 
both sides. Thus, it can not obey the periodic boundary 
condition with the same eigenvalue as it had for h < l/£k- 
its eigenvalue has to change! In practice, when we increase 
h the eigenvalue uj 2 starts to change rapidly and collide 
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Fig. 5: (a) DOS of our packings for 6 different pressures con- 
firming the main features of earlier studies [3-6,9] close to jam- 
ming, (b) Our frequency binned and angularly averaged values 
£(u))/L are all very similar. (c,d) Level spacing statistics for 
the modes that have £ > L in (c) and for the modes with 
£ < L. The lower frequency modes are essentially all extended 
and do show level repulsion in accord with the predictions from 
RMT [16,18], while the high frequency modes are truly local- 
ized and their level spacing is close to Poissonian. The gray 
lines indicate the frequency ranges used to obtain the level 
statistics in (c) and (d). 

with a neighboring eigenvalue when h w beyond that, 
when h > the eigenvalue w\ moves into the complex 
plane [17]. Hence we can simply obtain the localization 
length £fc of each mode k from the value hk at which the 
eigenvalue moves into the complex plane upon increasing 
h: £fc = 1/hk- Note that in this method we do not need 
to calculate the eigenfunctions explicitly — we only need 
to track the eigenvalues! 

It is straightforward to extend this method to higher di- 
mensions: as above, we simply multiply the off-diagonal el- 
ements of our dynamical matrix with an exponential e ry , 
where fij is the vector pointing from the center of particle 
i to its neighbor j. Our probe field h is now a vector, so by 
changing the angle that h makes with the x-axis, we can 
extract the angular anisotropy of the localization length 
£,((/)) of each mode. 

Results. — We first discuss some properties of the 
localization length of individual modes before turning to 
their scaling as a function of frequency, system size and 
distance from the jamming point. 

Anisotropy — Fig. 3 shows the angular dependence 
of three typical modes. One clearly sees that £((/)) is a ir- 
periodic function and that the angular variation of £ (4>) is 
significant. While few modes, like the second one in Fig. 3, 
have a quadrupolar structure, the anisotropy is predomi- 
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Fig. 6: Scatter plot for the localization lengths £ (determined 
to a precision of order unity) on a logarithmic frequency scale 
at p = 3 • 1CF 2 in (a) and p = 1CF 6 in (b) at system size 
L = 45. Note that at the large pressure the lowest-frequency 
modes are localized; at small pressures this is not the case. 
The inset illustrates the two lowest-frequency modes, which 
has (/Lps 0.3 in (a) and £/L « 2 in (b). 

nantly dipolar, as the histogram in Fig. 4(a) shows. 

We will denote from here on the angularly averaged 
value of the localization length of an individual mode by 
£. Figure 4(b) shows that the root mean square average 
angular variation A£ of £(</>) is almost half £, and that it 
is slightly larger at higher frequencies. There is no strong 
dependence of the anisotropy on the pressure, i.e., on the 
distance from the jamming point. 

Spread — The angularly averaged values ^(tu) also show 
a large spread, as Fig. 2 illustrates for a small value of 
the pressure. One also sees from this figure that most 
modes have a value of £ > L, which means that they are 
extended within the systems we can analyze — only our 
largest frequency modes are truly localized [20,22]. 

DOS — We now turn to a more systematic analysis of 
our data as a function of pressure and system size. In 
Fig. 5(a) we show that the density of states (DOS) of our 
packings behaves as found before [3-6, 9] for such pack- 
ings: As the the jamming point is approached by lowering 
the pressure, the density of low-frequency modes increases 
dramatically, which, as mentioned before, is due to the 
nearness of the isostatic point . 

Average localization length £(w) — For each dataset of 
the individual angularly averaged values of £, as in Fig. 2, 
we determine the frequency binned average values S;(uj) 
(each based on about 100 to 200 modes). The behavior of 
£/L as a function of (scaled) frequency is show in Fig. 5(b) 
for six different values of the pressure. In these average 
values, there is no strong variation with pressure, i.e. with 
distance to jamming. 

We already noted in Fig. 2 that most of our eigenmodes 
have £ > L, i.e. are extended in our finite system. This is 
also clear from Fig. 5(b): at all but the largest frequencies 
we have, £ > L. There are indeed roughly three regimes 
present in Fig. 5(b). From high frequencies towards low 
frequencies, we first have a range of high-frequency local- 
ized modes, for which £ < L. These modes are always 
present at any pressure and are the high-frequency modes 
in which only a few (light) particles oscillate more or less in 
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Fig. 7: (a) Finite size scaling for p — 4 • 10~ and linear system 
size L ranging from 15 to 45, confirming that the extended 
states, where £ > L, scale with L, while the high-frequency 
modes are, within the statistical error, //-independent, (b) 
Scaling collapse of £ according to (2) for hexagonal lattices with 
identical springs but varying masses mi oc 7?f, as for spheres. 
The distribution of radii Ri is taken to be flat, and W is taken 
to be the width of the distribution in percent. 

anti-phase as in an optical mode (such type of modes gen- 
erally arise immediately when disorder is introduced into 
an ordered system). For intermediate-range frequencies 
there is a plateau in £. Finally for the lowest frequencies 
(in the frequency range where actually the excess modes 
appear in the DOS in Fig. 5(a) at low pressures), there is 
an indication of an upswing in £ for small tu. We find this 
upswing at low frequencies in all our data, also on perco- 
lation lattices [19], where it is even more pronounced. 

Quasi-localized low-frequency modes at high pressure — 
From the above data for the bin-averaged £, it would ap- 
pear at first sight that we see no signature of the near- 
ness of the jamming point. This, however, is not true: in 
Fig. 5 we show data obtained by averaging over 100-200 
modes. However, this averaging washes out systematic 
trends visible for the lowest frequency eigenmodes discov- 
ered by Vitclli, Xu et al. [20,22]. When plotted on a 
logarithmic scale, as in Fig. 6, we see a systematic trend 
for £ of the low frequency modes to decrease with increas- 
ing pressure. As the inset of Fig. 6(a) illustrates, these are 
"quasi-localized" modes in which a reasonably well defined 
"localized" group of particles performs what looks like a 
resonant oscillation that is weakly coupled to the extended 
elastic field. For our limited range of L, we find £/L ~ 0.3 
and a reduced anisotropy of A£/£ m 0.2 for these modes. 

As we discussed in the introduction, for packings closer 
to the jamming point (at lower pressures) the isostatic- 
ity length £* increases as 1/AZ, where AZ is the excess 
contact number. Up to this scale I* we do not expect local- 
ized modes at low frequencies, since up to this scale the re- 
sponse mirrors that of the global floppy modes that emerge 
near the isostatic point. Indeed, within the system sizes 
we can study there are no low-frequency "quasi-localized" 
modes at all at low pressures, as Fig. 6(b) illustrates for 
p = 10 -6 , even though the response is in many ways more 
disordered due to the nearness of the jamming point! 

While our data are qualitatively in accord with the 
above scenario, we have unfortunately too few low- 
frequency "quasi-localized" modes to confirm quantita- 



tively that as we tune the packings closer to jamming the 
extent of the resonant region increases with £* ~ 1/AZ. 

Level Spacing Statistics — Based on the results of RMT 
[16, 18], one expects the following: the frequencies u>i of 
the localized modes should be independently distributed, 
so that their spacing Au>i = Wj+i — u>i obeys a Poisson 
distribution, while the modes which extend throughout 
the system should interact and repel each other, with a 
level spacing distribution which is given by the Wigncr 
surmise, Pw(s) = 7rs/2cxp (— 7rs 2 /4), where s = Aw/Acj. 
Figs. 5(c,d) confirm that this expectation is fully born out 
by our data at all pressures. Note that the distribution 
in Fig 5(c) deviates somewhat from the Wigncr surmise 
at the two highest pressures — this is due to the "quasi- 
localized" low-frequency modes discussed above. 

Scaling with system size L — One would of course ex- 
pect £ for the modes which extend through our system size 
to scale with L. As we will sketch below, we have used 
RMT [16,18] to derive for our method the scaling behavior 
£ ~ L d / 2 . More generally we propose 

I ~ L d/2 /W, (2) 

where W is a measure of the effective disorder. Fig. 7(a) 
shows that the £ ~ L-scaling is well obeyed for our two- 
dimensional packings for the extended modes in the range 
u> < 3 (as noted before the quasi-localized modes obey £ ~ 
0.3L), while the high-frequency localized modes for u> > 
3.4 have £'s which are indeed essentially L-independent. 

For our gently prepared granular packings the strength 
of the disorder can not easily be varied. In order to test our 
scaling prediction (2), we have prepared ordered hexago- 
nal lattices with all spring constants the same but varying 
masses. As Fig. 7(b) shows, we obtain very good data col- 
lapse with (2) at all but the highest frequencies. Note also 
that for small amount of disorder, we have £ 3> L. Re- 
sults for one-dimensional chains arc fully consistent with 
the predicted L 1 / 2 scaling [19]. 

Implications from Random Matrix Theory. 

Let us finally summarize what Random Matrix Theory 
(RMT) brings to bear on the study of the localization 
length. We refer for a more extensive discussion to [19]. 

1 ) In RMT it is well known that for analyzing the level 
statistics, like in Fig. 5(c,d), it is important to use the 
proper variable. The procedure to obtain the proper vari- 
able is the so-called "unfolding of the spectrum" [16]. For 
our case, the unfolding ensures that in each small interval, 
the mean level spacing is the same as in the original spec- 
trum. The proper variables are then indeed the frequen- 
cies ui, not the eigenvalues (J 1 of the dynamical matrix. 

2) The scaling (2) of the modes with £ > L can be 
understood as follows. When we turn on h, energy levels 
start to move on the real axis, some getting closer together 
and some further apart. Because of reflection symmetry 
under h — > —h (which is also apparent in Fig. 3), the shift 
is quadratic in h. We determine £ of a mode from the 
collision value h c at which two modes collide along the 
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real w-axis and "pop" into the complex plane. Accord- 
ing to RMT [16], the typical collision parameter is then 

h 2 mean level spaing p or systemS the mean level 

c typical level velocity J 

spacing is proportional to 1 /L d and the typical level veloc- 
ity does not depend on L d , from which the scaling £ ~ L d l 2 
immediately results upon identifying h c with 

3) In line with the large spread in the values of £, we 
find that the distribution P(£/£) implied by Fig. 2 falls 
off as (£/£) -3 f° r large £ both in 1 and 2 dimensions. 
This power law decay can be derived from how the level 
repulsion of the extended modes changes, when we change 
the perturbation parameter h by Ah [19]. 

Conclusions and Outlook. In this paper we have 
introduced a new method, motivated by previous studies 
of non-Hermitian quantum problems [17], which allows an 
analysis of localization in phonon spectrum, including the 
regime £ > L when the eigenmodes are extended within 
the finite systems we can study. The method is especially 
relevant for granular packings, where £ > L throughout 
most of the frequency range, since even in this regime our 
method gives different results depending on the amount 
of disorder. For the system sizes that are numerically ac- 
cessible at present we can not yet test the proposed scal- 
ing relation £ > I* ~ l/AZ quantitatively. Nevertheless, 
the disappearance of the "quasi-localized" low-frequency 
modes as we approach the jamming point by lowering the 
pressure, agrees with the scenario advanced by Wyart et 
al. [4-6] that up to this length scale the low-frequency 
rearrangements and modes extend over a diverging scale 
£*. We aim to study larger systems and more packings in 
the future using sparse matrix eigenvalue routines rather 
than a Mathematica program, to investigate the nature 
and scaling of the low-frequency modes in more detail. 

A few final remarks are in order, (i) Our method will al- 
low us to determine which type of disorder (mass disorder, 
bond disorder or topological disorder) plays the dominant 
effect in the localization behavior, (ii) The resonant region 
of the quasi-localized mode shown in Fig. 6 has quadrupo- 
lar symmetry, and is reminiscent of the quadrupolar de- 
formation fields that have been proposed [23] to dominate 
quasistatic shear relaxation. Although this is not true for 
all quasi-localized modes, this may not be accidental. The 
possible connection between these resonances and shear 
transformation zones is extremely intriguing and should 
be pursued further. (Hi) The results for in finite sys- 
tem typically show an upswing for small to, except at the 
largest pressures; whether this is a finite system analogue 
of the well known u> — > divergence in infinite 2d sys- 
tems [13] is unclear to us. (iv) The states with large but 
finite localization lengths at low frequency that we find 
at high pressures (see Fig. 6(a)) are intriguing. It will be 
interesting to see if these states persist in the presence of 
the entropic interactions at finite temperature, (v) Dif- 
fusion on percolation lattices is also an appealing model 
system to apply the method to: close to the percolation 
threshold most eigenmodes are truly localized and thus 



have £ <C L, while away from the percolation threshold 
there is a crossover to the regime where £ > L [19]. 
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